\chapter{微分方程}
\section{导论}
\[ a_1\cdot \dot{y(t)}+a_2\cdot y(t)+x(t)=0 \]


如果$ x(t) $是个常数，则该方程是自治的;如果$ x(t)=0 $,则该方程是齐次的。
解该方程有三种方法：
\begin{itemize}
	\item 图解法，可用于线性和非线性，但只适用于自治方程。
	\item 分析法，只适用于线性函数或是特定函数。
	\item 泰特展开，适用于非线性，但只能得到邻域内的情况变化。
\end{itemize}

\section{单方程}
\subsection{图解法}

考虑一个自治方程
\[\dot{y(t)}=f[y(t)]  \]

或者将其表达为两个不同的函数之和，
\[ \dot{y(t)}=f[y(t)]-g[y(t)] \]

首先找到稳态，即$ \dot{y}=0 $的点，然后，观察当$ f $位于$ g $上方时，$ y $要增大，$ f $位于$ g $下方时，$ y $要减小。该图纵轴是$ \dot{y} $,横轴是$ y $。

\subsection{分析解}
一般使用Matlab就可以解了。也可以用代数的方法求解，对于一个一般的常微分方程，
\[ \dot{y}(t) + ay(t)+ x(t)=0 \]
也即，
\[ \dot{y}(t) + ay(t)=- x(t)\]
\begin{itemize}
	\item 两边同乘以$ e^{at} $(这会使得下式左边很容易被积出来)，得到，
	\[ e^{at}[\dot{y}(t) + ay(t)]=- e^{at}x(t)\]
	\item 两边积分，得到，
	\begin{align*}
	e^{at}y(t)&=-\int e^{at}x(t)dt + C\\ 
	\Longrightarrow\hspace{2em}	  y(t)&=-e^{-at}\int e^{at}x(t)dt + e^{-at}C
	\end{align*}
	其中，$ C $为任意常数。
	
\end{itemize}

\section{线性常微分方程组}
一般可以书写如下，
\begin{align*}
\dot{y_1(t)} & =a_{11}y_1(t)+\cdots+a_{1n}y_{1n}(t)+x_1(t)\\
& \cdots\\
\dot{y_n(t)} & =a_{n1}y_1(t)+\cdots+a_{nn}y_{1n}(t)+x_n(t)
\end{align*}

上述方程组或者用矩阵形式写为，

\[ \dot{\bm{y}(t)} = \bm{A}\cdot \bm{y}(t)+\bm{x}(t)\]

\subsection{相位图}

该方法的缺点在于只适用于\underline{$2\times 2 $且具有稳态的自治的系统} 。相位图的横轴和纵轴分别为$ y_1,y_2 $，相位图中的箭头表示随着时间的流逝，变量的走向。

一个\textbf{对角线系统}的方程组如下，
\begin{align*}
\dot{y_1(t)} & =a_{11}y_1(t)\\
\dot{y_2(t)} & =a_{22}y_2(t)
\end{align*}

很容看到，
\begin{itemize}
	\item 系数全部为正，系统是不稳定的；因为其解非常直观，为，$ y_1(t)=e^{a_{11}x}+C_1,y_2(t)=e^{a_{22}x}+C_2 $。
	\item 全部为负，系统是稳定的；
	\item 一正一负，系统是鞍形稳定的
\end{itemize}

鞍形稳定存在一条稳定路径（稳定臂）和一条不稳定路径（非稳定臂），两者都经过稳态，实际上这个臂也是动态规划里面经常提到的政策函数。

主要来看看\textbf{非对角线系统},如下，
\begin{align*}
\dot{y_1(t)} & =0.06y_1(t)-y_2(t)+1.4\\
\dot{y_2(t)} & =-0.004y_1(t)+0.04
\end{align*}

比较容易地可以通过Matlab的$quiver(y_1,y_2,\dot{y_1},\dot{y_2})$函数来画向量场。或者按照传统的方式：第一步，画出$ \dot{y_1}=0,\dot{y_2}=0$两根线；第二步，考察两根线两边的情况；第三步，可以通过特征向量画出稳定臂和非稳定臂。具体的，看图\ref{sed}，教你一分钟画出相位图：
\begin{figure}[H]
	\includegraphics[scale=0.37]{aa.png}
	\caption{怎么画相位图}\label{sed}
\end{figure}


实际上，非对角线系统可以通过对角化转化为对角线系统，那么对应的，
\begin{itemize}
	\item 两个特征值全部为正，系统是不稳定的；
	\item 两个特征值为负，系统是稳定的；
	\item 两个特征值一正一负，系统是鞍形稳定的
	\item 两个特征值是实部为负的复数。以震荡方式向稳态收敛。
\end{itemize}

为何复数是这样？参见时间序列分析中差分方程的解法。首先意识到任何高阶微分方程可以通过变量替换而改写成一阶微分方程组。因此，对于双变量的一阶微分方程组，可以考察一个二阶微分方程。比如，对于一个如下的二阶常系数微分方程，
\[ y''+py'+qy=0 \]

尝试$ y=e^{rx} $是上述方程的解。将这个猜测的解代入上式，因为$ y''=r^2e^{rx},y'=re^{rx} $，所以有，
\begin{align*}
(r^2+pr+q)e^{rx}=&0\\
\Longrightarrow r^2+pr+q=&0
\end{align*}

因此，只要$ r $满足上式，那$ y=e^{rx} $就是这个方程的解。上式通常称为微分方程的特征方程。特征方程的两个根有三种不同情形，不等的实根，相等的实根以及共轭复根。因为，
\[ r_1=\frac{-p+\sqrt{p^2-4q}}{2},r_2=\frac{-p-\sqrt{p^2-4q}}{2} \]

先来考察共轭复根的情况，即
\[  r_1=\alpha+\beta i,\hspace{2em}r_2=\alpha-\beta i, \]
其中$ \alpha=-\frac{p}{2},\beta = \frac{\sqrt{4q-p^2}}{2}$。此时方程的解为$ y_1=e^{(\alpha+\beta i)x},y_2=e^{(\alpha-\beta i)x} $。为消掉那个$ i $，利用欧拉公式$ e^{i\theta}=\cos \theta+i\sin\theta $，则，
\begin{align*}
y_1=&e^{(\alpha+\beta i)x}=e^{\alpha x}e^{i\cdot \beta x}=e^{\alpha x}(\cos\beta x+i\sin\beta x)\\
y_2=&e^{(\alpha-\beta i)x}=e^{\alpha x}e^{-i\cdot \beta x}=e^{\alpha x}(\cos\beta x-i\sin\beta x)
\end{align*}

因为对于一个二阶微分方程，如果$ y_1,y_2 $是它的两个解，那么$ y=C_1y_1+C_2y_2 $也是它的解。所以把上述两个复值函数相加除以2得到的实部和相减除以$ 2i $得到的虚部，也可以看作"$ y_1,y_2 $"，从而得到，
\begin{align*}
\bar{y}_1=\frac{1}{2}(y_1+y_2)=e^{\alpha x}\cos \beta x\\
\bar{y}_2=\frac{1}{2i}(y_1-y_2)=e^{\alpha x}\sin \beta x
\end{align*}
它也是方程的解。注意到此时方程的解就呈现出周期性。如果实部为正，则震荡发散；实部为负，则震荡收敛；实部为0，则就是一个椭圆($ y_1,y_2 $分别为横、纵坐标)。

\subsection{分析法}

无论是齐次（$ x(t)=0 $）还是非齐次($ x(t)\ne 0 $)，该法的一大特征在于将方程组中的每个方程的右边转化成只包含自身的函数,这样每个方程又可以应用单方程的分析法来求解。如对于齐次的而言（仅为表述方便），就是把$ A $对角化，即找一个方阵$ V $，使得，
\[ V^{-1}AV=D \]

然后有
\begin{align*}
\dot{z(t)} & \equiv V^{-1}\dot{y(t)}\\
& =V^{-1}Ay(t)\\
& = V^{-1}AVV^{-1}y(t)\\
& =D\cdot z(t)
\end{align*}

余下求解显而易见。


\section{非线性方程组}
\subsection{泰勒展开}

通过将方程组在稳态附近线性化，则可以按照上面的方法继续求解。
\subsection{时间消去法}
对于一个如下的非线性方程组，
\begin{align*}
\dot{k} & =k^{0.3}-c\\
\dot{c} & =c(0.3k^{-0.7}-0.06)
\end{align*}

其边界条件为$ k(0)=1,lim_{t\rightarrow\infty}[e^{-0.06t}k(t)]=0 $。如果我们还知道$ c $的初值，那么通过上述方程组不断迭代，就可以得到整个$ k $和$ c $的路径。很遗憾，我们只有$ k $的初值。但要得到策略函数还是有技巧的，即通过如下技巧，
\[ c'(k)=\frac{\dot{c}}{\dot{k}}=\frac{c(0.3k^{-0.7}-0.06)}{k^{0.3}-c} \]

这是一个关于$ k $求导的常微分线性方程，而且可以把稳态值作为初值进行求解。如果是从数值上进行追踪的话，困难在于作为初值的稳态点$ c'(k)=\dot{c}^*/\dot{k}^*=0/0 $,即一开始就迭代不了了，该问题要么通过洛必达法则解决，要么意识到稳态时策略函数的斜率等于负特征向量的斜率。然后从该点开始，逐步得到整个$ c(k) $。
\section{最后的一点体会}
上面的这些解法，通常对于小规模比如2个状态变量的玩具模型是非常合适的。模型一旦复杂化，即变量数目增加到3个以后，通常的做法是对非线性微分方程组进行对数线性化，然后再来求解。此时的求解方法，通常是数值的，可以参考《RBCofABC》。